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Abstract 

We present in detail a technique to extract the potential between a 
static quark and anti-quark pair from Wilson loops measured on dynamical 
configurations. This technique is based on HYP smearing and leads to an 
exponential improvement of the noise-to-signal ratio of Wilson loops. We 
explain why the correct continuum potential is obtained and show numerical 
evidence that the cut-off effects are small. We present precise results for the 
non-perturbative potential. As applications, we determine the scale r^ja 
and study the shape of the static potential in the range of distances around 
ro, where it can be compared with phenomenological potential models. 
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1 Introduction 



The potential V (r) between a static (infinitely massive) quark and anti-quark pair 
separated by distance r can be computed from lattice quantum chromodynamics 
(QCD). It is extracted from the expectation values of Wilson loops, which are 
traces of products of links along rectangular paths extending in Euclidean time 
and one spatial direction. In this article we consider only on-axis Wilson loops 
but off-axis (non-planar) Wilson loops can also be used. Alternatively the static 
potential can be extracted from the correlator of two Polyakov loops. Due to 
confinement, the signal of Wilson loops falls exponentially with the area of the 
loop (until string breaking sets in) but their variance is approximately constant. 
In the statistical average of standard Monte Carlo lattice simulations, the signal 
of Wilson loops is the result of strong cancellations between positive and negative 
contributions. This leads to an exponentially growing noise-to-signal ratio which 
prevents the calculation of the potential at large distances. 

In pure gauge theory this problem has a cure. An exponential suppression of 
the statistical noise of Wilson loops can be achieved by the multi-hit (or one-link) 
method [1] and much further by the multilevel algorithm [2J. These algorithms 
are not applicable in presence of dynamical fermions due to the non-locality of the 
effective gauge action when the logarithm of the fermion determinant is included. 
In [5] a smearing technique called hypercubic (HYP) smearing was introduced 
which can also be used to measure Wilson loops in the presence of dynamical 
fermions [3] • In pure gauge theory it was demonstrated in [5] that the use of HYP 
smeared links leads to a determination of the static potential comparable in preci- 
sion to the multi-hit method. In [6] a new action for static quarks was proposed, 
which uses HYP smeared links in the time covariant derivative of the Eichten- 
Hill action. This leads to an exponential reduction (compared with using the 
Eichten-Hill action) of the noise-to-signal ratio for heavy-light correlation func- 
tions. This effect is due to the fact that HYP smearing in the static action reduces 
the coefficient of the divergent part of the self-energy of a static quark 

The interest in the determination of the static potential V(r) through lattice 
simulations is twofold. On the one side, there is the possibility to set the scale 
(i.e., determine the lattice spacing) through the scale r introduced in [9]. The 
latter is defined from the static force F(r) = V'(r) as the solution of 



A physical value for the scale r$ (0.45 . . . 0.5) fm can only be determined through 
phenomenological potential models. It is desirable for an absolute determination 
of the lattice spacing to use a quantity which is directly accessible from experiment 
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and replace tq by a quantity like a hadron mass or decay constant. But still, ro is 
very useful for a relative scale setting. 

On the other side, the static potential is an interesting observable by itself 
for phenomenology (see the conclusions) and to study the structure of gauge 
theories [T0l - fl3] . It exhibits clear effects of dynamical fermions, such as string 
breaking [14] . see the latest study in QCD [15] and high precision studies with 
multilevel algorithms in other models p.6j. In order to study the potential at the 
distances where the string breaks, operators which dominantly create static-light 
meson pairs have to be included in addition to the Wilson loops and we plan to 
do so in the future. In this article we will concentrate on the determination of the 
static potential at distances smaller than the string breaking distance r b « 3r 
[T7] . We will study the quantity 

c(r) = ^r 3 F'(r). (1.2) 

It is a physical, renormalized quantity, which can be used to define a running 
coupling. In [TS] c(r) has been determined with high precision in pure gauge 
theory using a multilevel technique. We will compute it in this article for QCD 
with N{ = 2 flavors of quarks. 

In section two we will describe our techniques to extract the static potential 
from HYP smeared Wilson loops. We explain why this procedure leads to a de- 
termination of the continuum static potential up to 0(a 2 ) lattice artifacts, which 
appear to be small. In section three we present our results for the static poten- 
tial, the scale r^/a and the quantity c(r) determined on a configuration ensemble 
generated with Wilson gauge action and iVf = 2 flavors of 0(a) improved Wilson 
quarks at /3 — 5.3. The quark mass corresponds to a pseudoscalar mass value 
close to r m PS = 1 and we get a value r /a = 6.75(6). 

2 Techniques 

2.1 Static potential with HYP smearing 

We measure r/a x T/a on-axis Wilson loops W(r,T) on gauge configurations 
generated with Nf = 2 dynamical fermions. The technique is based on HYP 
smearing and was introduced in [3]. Before measuring the Wilson loops, we replace 
all the gauge links by HYP-smeared ones. We consider two choices of the HYP- 
smearing parameters: one is 

oti = 0.75, a 2 = 0.6, a 3 = 0.3, (2.1) 
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Figure 1: Schematic representation of the measurement of Wilson loops. In the first 
step (left figure) only the spatial Wilson lines are HYP-smeared: this corresponds to 
the definition of an operator that creates a \QQ(r)} state. In the second step (right 
figure) the temporal Wilson lines are HYP-smeared: this corresponds to the choice of 
the static quark action (and a modification of the operator O). 

which we refer to as HYP, and the other is 

ai = 1.0, a 2 = l-0, a 3 = 0.5, (2.2) 

which we refer to as HYP2. We adopt the approximate projection onto SU(3) as 
described in [7] and always use Eq. (2.24) and four iterations of Eq. (2.25) in [7J. 

In the following we show that this procedure leads to a determination of the 
potential between quark and anti-quark sources that agrees with the continuum 
potential up to 0(a 2 ) effects (after renormalization). The ingredients in this 
demonstration are the selfadjoint positive transfer matrix of the lattice gauge 
theory with Wilson fermions and Wilson plaquette action (rigorously proven [T9] ) 
as well as the existence and universality of the continuum limit of the lattice 
theory with a static quark (lowest order of heavy quark effective theory [20]). The 
latter property is generally assumed and has been tested frequently (see [21] for 
a longer discussion). 

For the purpose of showing Eq. fl2.6p . we split the HYP-smearing of the links 
used in building a Wilson loop into two steps, which are schematically represented 
in Fig. [TJ In the first step we consider Wilson loops where only the space-like links 
are HYP smeared. The smearing involves links at time-slices^ t = —a, 0, a and 
t = T — a,T,T + a and corresponds in the Hamiltonian formalism to an operator 
and O that creates or annihilates a state \ip®Q(r)) consisting of a static quark 

1 In total there are N t time-slices and periodic boundary conditions are imposed in all 
directions. 



and anti-quark pair at time-slices t = a and t = T — a respectively. The static 
sources are separated by a distance r along one of the spatial directions. The 
path integral average of this Wilson loop can be written as a quantum mechanical 
expectation value 

Tr \T^- T/a - 2 OT q g(r) T / a - 2 d^ } 
(W(r, T)) = ^ } - , (2.3) 

where To is the transfer matrix in the gauge-invariant (or zero charge) sector 
of the Hilbert space, Y qq -{r) the transfer matrix in the sector with a quark and 
an anti-quark source separated by r and Tr is the operator trace in the Hilbert 
space. We denote the transfer matrix in the temporal gauge (where the time- 
like links are set to identity) by T tcmp . The Hamiltonian HI is defined through 
aM = — ln{T temp }. For the theory with Wilson quarks without a clover term 
Ttemp has been constructed in [19] . The transfer matrix in a specific charge sector 
is obtained by multiplying T temp with the projectors onto that charge sector. 
Note that the representation Eq. (12. 3p differs from the usual one only in that the 
operators O represent fields in the path integral on three timeslices, not one. If 
written down explicitly in the form of [19] they involve integration kernels. But 
their explicit form is not relevant here. Using the spectral decomposition of the 
transfer matrices (see for example [22J) and taking the limit N t — > oo, Eq. (12. 3 p 
becomes 

(W(r,T)) Nt ^ Y, c ^>~ Vn(T)(T ~ 2a] '» ( 2 - 4 ) 

n 

where c n = (n; qq\O^\0) are the overlaps of states obtained by applying the oper- 
ator O' to the vacuum |0) with the eigenstates \n;qq) of the Hamiltonian (with 
eigenvalues E n (r)) in the sector with a quark and an anti-quark source. In 
Eq. ( 12. 4 p we use V n (r) = En(r) — E^\ where E^ is the vacuum energy. For 
example V (r) is the static potential and V\(r) its first excitation. 

In the second step we rewrite the Wilson loop as a path integral expectation 
value 

(W(r,T)) = ~(^ h {a,6)P4a,0;a,rk)^- h {a,rk) 

$- h (T - a,rk)Pl(T - a, 0; T - a,rk)^ h ( T ~ M)) , (2.5) 

where i^h^jph an d ipRi are ^ ne static quark and anti-quark fermion fields re- 
spectivel)o and P±(t, 0; t, rk) represents the gauge parallel transporter made from 

2 The prefactor — \ and the gamma- matrices are due to our choice of treating the static 
quark fields as 2-component static fermion fields, see for example |21j . 
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a product of space-like HYP-links at time t ± a and temporal links at time t ± a 
in =F-direction (dashed lines in Fig. CQ. 

After integration over the static fields, the static quark propagator generates 
the time-like links in the observable, cf. Eq. (2.4) in [7], and one recovers the Wil- 
son loops. Different choices for the static quark action can be made, in particular 
we consider here the one where the covariant derivative in time involves HYP or 
HYP2 linkaH [6J. It was shown in [23] that static potential differences (where the 
self energy is canceled) have 0(a 2 ) leading lattice artifacts, essentially due to the 
automatic 0(a) improvement of heavy quark effective theory [24J. This is true in 
the theory with dynamical fermions provided that they are 0(a) improved^ We 
therefore conclude 

K HYP/HYP2 (r) _ 2 £HYP/HYP2 = ^continuum (f>) _ 2 £-nti_ + 0(fl2) ? q g) 

where for convenience we have subtracted V(oo) = 2E stat . Here E stat is the 
binding energy of a meson made of a static and a light dynamical quark. 

In order to investigate the magnitude of the lattice artifacts we compare in the 
left panel of Fig. [2] the qq-coupling a qq (r) = r 2 F(r)/Cp for the HYP and HYP2 
actions. The static force F(r) is obtained from the static potential F(r — a/2) = 
[V(r) — V(r — a)}/ a. Details about the extraction of the static potential from 
correlation functions of Wilson loops are presented in the next section. We use the 
dynamical configurations described in Section |3j The difference of the couplings 
is given, to leading order in the cut-off effects, by a 2 /r 2 G(Ar,mr). The function 
G describes the r-dependence and quark-mass-dependence of the cut-off effects. 
The size of the cut-off effects is small but with our errors they are significantly 
different from zero for r < r . They happen to be most significant at r w r /2. 

In Section I2TB1 we describe how improved observables can be defined such that 
these cut-off effects are eliminated at tree level and are substantially reduced non- 
perturbatively [23]. The right panel of Fig. [2] is the same as the left but using the 
improved definition of the force Eq. f)2.14p . The cut-off effects are visibly reduced. 
We emphasize that the figure is not sufficient to exclude cutoff effects which are 
independent of the choice of static action. Different lattice spacings are needed 
to study those. 

Our choice of the static quark action and the smearing of the spatial links 
is with parameters HYP2 Eq. (12. 2p . It gives a static potential with a somewhat 

3 We smear also the temporal links contained in the definition of the parallel transporters 
P± in Eq. (|2.5p . This corresponds to a change in the definition of the operator O in Eq. (I2.3|) 
and has no consequences for the argument presented here. 

4 For Wilson fermions improvement is achieved by adding the clover term [25H27] or by using 
a twisted mass term [25J "at maximal twist" [29] , 



6 









■ 




- 


+ 

e 






■ 9 


+ 
o 


HYP 
HYP 2 


0.4 0.5 


0.6 0.7 


0.8 



r/r 




Figure 2: The qq-coupling a qq (r = 1///) obtained from the static force F(r — a/2) = 
[V(r) — V(r — a)]/ a with two different choices of the static action (left panel). In 
Section 12.31 we give an improved definition of the force which is free of cutoff effects at 
tree level of perturbation theory (right panel). 



better statistical precision than with parameters HYP. This can be understood 
in perturbation theory: the HYP2 parameters are such that they approximately 
minimize the one-loop coefficient of the 1/a self-energy contribution of a static 
quark [318]. Our data show that this property remains true non-perturbatively: 
we find V KYP2 - V nYP « -0.07/a « 2(£ S ^ P2 - E™). For the last statement 
we use the results for P stat of reference [30] . 

2.2 Variational basis 

On the HYP2-smeared gauge link configurations {U(x,fi)}, we measure a corre- 
lation matrix of on-axis Wilson loops at fixed spatial extension r/a and temporal 
extension T/a: 

C lm (T) = (tr{p (/) (0;rA;)P(rA;;rA; + T6)P (m)t (r6;rA; + T6)P t (0,r6)}\ , 

(2-7) 

where P(x,y) represents the product of links connecting y to x. Neglecting the 
superscripts on the spatial P's, Eq. (12. 7p is equivalent to Eq. (I2.5P after integrat- 
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ing out the static fields. In the product of spatial links, the superscript means 
that the links Ui(x,k) used in the product are obtained by applying the spatial 
smearing 5 s hyp operator rii times 

Ufak) = (S s}lYP pU(x,k). (2.8) 

5 s hyp means smearing with only two levels of HYP blocking with staples re- 
stricted to spatial directions and therefore it needs two parameters, which we 
set to «2 = 0.6 and a 3 = 0.3. In the argument of the previous section the fat 
parallel transporters P"' correspond to operators 0[ implementing trial states 
\"4>i {r)) = Oj\0). In [3T] a formula for suitable smearing parameters rii is given 
in the case of APE smearing. In order to choose n\ for spatial HYP smearing, 
we use the result of [32], that the mean squared extension of APE smearing is 
approximately a ape ^z,ape a 2 /3, and require that this is equal to n/a 2 for HYP 
smearing. We get an approximate formula for a good range of HYP smearing 
levels 

n * (2.9) 



12 V a J 

For our data on the configuration ensemble E5g (see Section [3]) we have computed 
a large correlation matrix using smearing levels no,i, 2,3,5- We find that this basis 
can be reduced to an optimal subset of M = 3 levels 

n 2 = 8, 713 = 12, n 5 = 20, at = 5.3. (2.10) 

The higher smearing levels improve the determination of the energy levels. 

We use the generalized eigenvalue method [3T1I33TI35] to extract the ground 
state potential as follows. We first solve the generalized eigenvalue problem 

C(t)4> a = \ a (t,t )C(to)ip a . (2.11) 

Then we perform a fit to 

E a (t + - * ) = In (X a (t, t )/X a (t + a, t )) = E a + ^-(^-^W+f ) , ( 2 .12) 

with fit parameters E a , (3 a and E M , simultaneously for a = 0, 1 (a = corre- 
sponds to the ground state, a = 1 to the first excited state), to/a = 2,3,4 and 
to + a < t < 2to (the latter constraint is necessary for eq.( 12.12l) to hold [35] ), i.e., 
we have 18 data points E a (t + |,t ) for 5 fit parameters. 

The values of the ground state potential V(r) as a function of r are determined 
from a plateau average of the corresponding effective masses Eo(t,to) starting at 
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Figure 3: Effective masses Eo(t,to) (filled blue circles for to /a = 5, empty black dia- 
monds for to/ a = 2) for the ground state potential at r = 7a. The red dotted line is the 
fit Eq. (|2.12p . The blue line is the plateau average from the points in the blue shaded 
area (the blue dashed-dotted lines are the plateau errors). 



the value t = 2t + |, where the fixed value t is determined by the requirement 
that 

cr sys (£ (2t + ~ i )) = /3 oe -(^-TO>+!) < ia stat (£ (2t + , (2.13) 

where cr sys (-) and <7 s tat(-) denote the systematic and statistical error respectively. 
For our data, Eq. (12. 13[) is satisfied for t /a = 5 for all values of r. The plateau 
average is stopped before the time, when either the difference of the effective mass 
with the one at t — 2t + | is larger than the statistical error of the latter or the 
statistical error of the effective mass is larger than twice the one of the effective 
mass at t = 2t + f • The effective masses E Q (t, t ) (filled blue circles) together with 
the fit Eq. (I2.12p (red dotted line) and the plateau average (blue line with error 
band marked by blue dashed-dotted lines) for r = 7a m r$ are shown in Fig. [3j 
The plateau average comprises three points at t/a = 10.5, 11.5, 12.5. The error of 
the plateau average is the sum of the statistical and the systematic errors, with the 
latter being given by the left-hand side of Eq. (I2.13p . For comparison, we also plot 
in Fig. [3] the effective masses obtained using to /a = 2 (empty black diamonds). 
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Figure 4: Effective masses E\{t,to) (filled blue circles for to/a = 4, empty black dia- 
monds for to /a = 2) for the first excited state potential at r = 7a. The dashed black line 
represents the value 2a£^ sta t j the meaning of the other curves is as explained in Fig. [3j 



They are part of the data set fitted using Eq. (I2.12p . At the times when they 
are both defined, the effective masses for to/ a = 5 and t /a = 2 agree with each 
other, which is somewhat surprising since non-leading corrections certainly have 
a dependence on to- 
la principle the excited potentials can be determined in the same way. How- 
ever, the analysis is complicated by the dynamics of string breaking. From model 
studies [3&H38] as well as from [15], we know that an extraction of the potentials 
requires the inclusion of operators which dominantly create static-light meson 
pairs in addition to the string-like operators we use here. Only then does the 
ground state at large distances r > contribute significantly to the spectral de- 
composition of the correlation function matrices at the accessible time separations 
(cf. [39]). While we are not concerned here with this string breaking region, it is 
known [THJGIE] that for r < the first excited state is an (approximate) meson- 
anti-meson state at V\ ~ 2E stat . This state is not well seen in our computation 
which does not include the meson pair operators. In Fig. H] we show the effective 
masses E\ for r = 7a ~ r . The dashed black line represents 2aE stat = 0.7007(14), 
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the meaning of the other curves is as explained for Fig. [3j Although the effective 
masses seem to form a plateau at times t = 8.5a and t = 9.5a, they drop at values 
t > 10.5a, but the statistical errors are too large in this region to determine an 
energy level. For r < r we see plateaus for E\ which are compatible with 2aE stat . 

This deficiency of our variational basis also affects the estimate of the ground 
state potential, but here the only concern is the description of the corrections to 
the asymptotic plateau of the form /3oe~^ EM ~ Eo ^. These enter the final numbers 
and errors only to estimate where we start the plateau average such that excited 
state contaminations are small compared to our statistical error. A very precise 
determination of Em or in general of the excited states is not necessary for this 
purpose. Furthermore, the effective mass figures show that our plateau selection 
is rather conservative; the extracted ground state potential is reliable within the 
cited errors. 



2.3 Tree level improved force 

In order to determine the scale r [3] from Eq. ( II. ip we will need the static force 
F(r). An improved definition of the force on the lattice is [9|l2"3"tl31~] 

F(n) = [V(r) -V(r -a)]/ a, (2.14) 

where rj = r — a/2 + 0(a 2 ) is chosen such that at tree level in perturbation 
theory [ID] one has 

F t ree(ri)=C F -^ 1 , (2.15) 
47T7j 

where Cp = 4/3 for gauge group 577(3). The formula for r\ depends on the static 
quark action and it is given in Appendix |A] for HYP actions. 

An improved lattice definition of the quantity c(r) in Eq. ( II ,2p is given by [IB] 

c(f) = ^ 3 [V{r + a) + V{r-a)-2V{r)}/a 2 , (2.16) 

where r = r + 0(a 2 ) is chosen such that 

c tree (r) = -C F |L (2.17) 

The formula for f depends on the static quark action and it is given in Appendix 
lAl for the HYP actions. In Appendix [B] we give the 4-loop beta function for the 
coupling a c = —c(r)/Cp which we will use to generate perturbative curves for 
c(r) to be compared with the lattice data. 
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Figure 5: Ground state potential V(r) (circles). The shaded area marks energy states 
larger than Em (here M = 3) determined from Eq. (|2.12p . 



3 Results 

We compute the static potential on the lattice ensemble E5g generated by the 
CLS ("Coordinated Lattice Simulations") project at /3 = 5.3, k = 0.13625 with 
geometry 64 x 32 3 and periodic boundary conditions for all fields apart from 
anti-periodic boundary conditions for the fermions in time. The value of the 
pseudo-scalar mass is amps = 0.15. The algorithm used in CLS is the deflation 
accelerated DD-HMC algorithm jJTJHS]. The trajectory length is r = 4 and we 
separate the measurements of Wilson loops by 4 trajectories. Given the block size 
8 4 , the active links represent 37% of all links. Hence the separation in molecular 
dynamics time between measurements is approximately 6 units (when all links are 
changed). We have a statistics of about 1000 measurements. 

In Fig. [5] we show the ground state potential V, and for illustration the 
rough estimate of the excitation Em (here M = 3) in Eq. (12.121) . In order to 
get renormalized quantities we subtract twice the binding energy E stat of a meson 
made of a static and a light dynamical quark. Everything is made dimensionless 

5 https : / / twiki . cern . ch/twiki/bin/view/CLS/WebHome 
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T\ja a 2 F{ri) r/a c(r) 



3.55805 0.05776(12) 

4.52674 0.04690(19) 

5.50073 0.04074(35) 

6.48362 0.03699(48) 

7.47397 0.03393(75) 

8.46922 0.0325(12) 

9.46734 0.0295(13) 

10.4670 0.0287(19) 

11.4676 0.0310(25) 

12.4685 0.0248(53) 

13.4697 0.0285(48) 

14.4709 0.0373(75) 

15.4721 0.030(11) 



4.046306 -0.3596(41) 

5.026094 -0.391(15) 

5.999703 -0.405(41) 

6.977869 -0.519(99) 

7.917429 -0.35(22) 



Table 1: The values of the force F(r) in lattice units and the physical quantity c(r) at 
the accessible improved distances n and r respectively. We do not include values that 
require the potential V(r = 2a), since it may be affected by large cut-off effects. 



by appropriate powers of ro extrapolated to the chiral limit, see below. The first 
excited state potential is not shown due to the difficulties described above. 

The range of string breaking is not yet reached. We can estimate it from the 
condition V(rb) = 2E sta .t to be 



2.4 < ^ 

— ro 



< 2.6. (3.1) 

ro mps=1.0 

For comparison, in [TJ] 7b/V ~ 2.5 was found at a larger quark mass corresponding 
to r mps = 1.7, albeit in the theory without 0(a) improvement. 

The scale ro is defined from the condition Eq. (11. ip . The static force is 
computed from Eq. (12.141) using the improved distance r\ in Eq. (lA.ip . In Table [1] 



we list the values of the force in lattice units. We do not include the force at 
ri/a = 2.58875 because it requires the potential at distance r = 2a, which may be 
affected by relatively large cut-off effects. We determine the solution of Eq. (11.11) 
by interpolation of the force F, using a 2-point interpolation F(r) = /o + /V?" 2 
and a 3-point interpolation adding a /4/r 4 term to control the systematic error 
(it is found to be negligible). We obtain 

= 6.747(59) (0 = 5.3). (3.2) 

amps =0.15 
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Figure 6: Auto-correlation function p(t) and integrated auto-correlation time T m t of tq. 
The Monte-Carlo time is in units of molecular dynamics time. 



The error is determined by taking the upper bound r int = 6 (see below) and 
neglecting the systematic error (due to excited state contributions), which is much 
less than the statistical one (due to condition Eq. (" 12 . 1 3[) ) . In jl3] we presented a 
preliminary value r$/a = 7.05(3) extrapolated to the chiral limit. Throughout this 
article we use this value (without errors) for the purpose of plotting dimensionless 
quantities. 

In jH] it was shown that the auto-correlation time of the topological charge 
suffers from critical slowing down proportional to a~ 5 in the present range of 
lattice spacings. However, in the same reference it was shown that Wilson loops 
are affected by a much milder critical slowing down, implying that their coupling 
to the slow modes of the Monte Carlo simulation is small. A method for correcting 
the error analysis, by adding a tail to the auto-correlation function that takes into 
account the coupling to slow modes, was presented in [44J . We use this method in 
our data analysis and we set r exp = 39 in molecular dynamics units^l from Table 
4 of [43]- In Fig. [6] we show the auto-correlation function p(t) and the integrated 

6 In these units, the DD-HMC molecular dynamics time is multiplied by the ratio of active 
links, which in our case is 37%. 
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Figure 7: The physical quantity c(r) in Eq. (jl.2|) . Comparison of iVf = 2 (circles) with 
Nf = (pluses) Monte Carlo data (for iVf = taken from [18J) and perturbation theory. 
Also the value c = —0.52 in the Cornell [46 1 potential and the curve (dotted) derived 
from the Richardson [47| potential are plotted. 

autocorrelation time Ti nt of ro, determined with the program@of [44] implementing 
the method of jHJHH]. The vertical dashed lines in the plots mark the applied 
summation windows, the lower one is used when we add the tail due to the slow 
modes, while the larger one comes from using the method of [45J. Adding to the 
summed autocorrelation function the correction due to the slow modes leads to 
the upper curve and upper bound on r int , which we take for all quantities in our 
analysis. The lower curve corresponds to r int determined from [45J. For r we 
get an upper bound r int = 6 which is a factor 6.5 smaller than r exp , but a factor 
1.5 larger than without accounting for effects of undetected slow modes (lower 
bound). 

In Fig. [7] we plot our result for the physical quantity c(r) computed from 
Eq. (I2.16P using the improved distance f in Eq. (1A.5j) . The numbers are given in 
Table[TJ In order to compare our Nf = 2 results (circles), we plot them together 
with the Nf = data (pluses) of [18] and with the perturbative curves obtained 

'http : / / www-zeuthen . desy . de/alpha/ 
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using the 4-loop beta function (continued line for Nf — 2, dashed-dotted line for 
N{ = 0). The perturbative formula for c(r) is presented in Appendix |B] and we 
used the preliminary updated value of the A parameter presented in [13]. The 
spread of the perturbative curve reflects the uncertainty of the A parameter. For 
a comparison with our Monte Carlo data, it is legitimate to plot the perturbative 
curve of c(r) in massless perturbation theory, since quark mass corrections are 
of order a 2 x (m q r) 2 and are expected to be negligible at our small quark mass. 
The distances in Fig. [7] are normalized by ro extrapolated to the chiral limit. As 
the perturbative curves already indicate, the value of c for Nf = 2 is found to 
be lower than for Nf = 0. In pure gauge theory, c(r) starts approaching the 
asymptotic value c(oo) = — n/12 with corrections of order 1/r as predicted from 
the effective bosonic string theory [48,49j. Our data for Nf = 2 have quite large 
errors when r/ro > 1. We compare them with the value c = —0.52 that it takes 
in the phenomenological Cornell potential |46j and with the curve obtained from 
the Richardson potential [57]. Our data seem to follow the Richardson curve for 
r < r quite closely. It is not yet possible to tell whether there is a plateau region 
around or above r before string breaking sets in. We will return to this quantity 
in our future studies. 

The comparison to the purely perturbative curve shows qualitative agree- 
ment. A meaningful quantitative comparison requires a careful study of lattice 
artifacts which may be quite noticeable in the region of small r, where perturba- 
tion theory applies. Indeed perturbation theory by itself suggests that at least 
r — f r o is necessary [50], in particular when the new 4-loop beta function is taken 
into account as discussed in appendix B. 

4 Conclusions 

We have presented a detailed analysis of the static potential defined by the HYP2 
action for the static quarks. Fig. [5] and Fig. [TJshow the quality of our data. Judged 
by a comparison of HYP and HYP2, cut-off effects in the potential appear to be 
small. The scale r^/a can be determined with precision better than 1%. We 
observe large effects due to dynamical fermions in the quantity c(r) defined in 
Eq. (JT2J. 

As can be seen in Table [IJthe error on the force grows faster with the distance 
r as compared to the pure gauge case (see Table 2 of [IS])- This effect is amplified 
by r 3 for the quantity c(r). It remains to be seen whether the inclusion of fermionic 
correlators in the variational basis will lead to an improvement due to a larger 
overlap with the ground state and the resulting earlier start of a plateau. 

A precise study of the static potential is relevant for phenomenology in an 
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indirect but important way. As reviewed in [51], there is an impressive effort 
to apply potential non-relativistic QCD (pNRQCD) [52] to the top - anti-top 
production in a future e + e~ collider and to many other processes. This effective 
theory includes ultrasoft gluons and is treated perturbatively in the QCD coupling. 
While the potential of pNRQCD is not the same as the static potential, the two 
are intimately related; they differ only starting at NNNLO accuracy. It is hence 
very useful to understand where the perturbative approximation to the static 
potential can be trusted. Fig. [7] is a start for that, but a precise investigation 
requires the removal of lattice artifacts [23]. In the future we plan to work both 
on this connection to the perturbative regime of QCD and on the large distance, 
string breaking, region. 
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A Improvement 



The tree level perturbative expression for the static potential, which is extracted 
from Wilson loops where the static quark line is HYP smeared, is given in [HI55]. 
From it we easily derive the formula for rj defined from Eq. (12. 15ft : 



2\-l 



with 



(47rrj 



-[G HY p(r, 0, 0) - G HY p(r - a, 0, 0)]/a . 



d 3 p U 3 j=l COS ( X jPj/ a ) X /sm(p) 



a (2vr)3 

where r — (xi, x 2 , x 3 ), pj = 2sin(j9j/2) and the smearing factor is 



fsmip) 



(A.l) 



(A.2) 



(A.3) 



a 2 



Q j0 (p) = l + a 2 (l + a 3 )- ^(l + 2a 3 )(pl+p 2 2 +Pl-p 2 j ) + -^ H Pi 



(A.4) 



17 



r/a n/a r/a 

HYP HYP 2 HYP HYP 2 



4 


3.48560 


3.55805 


3.97292 


4.04631 


5 


4.45369 


4.52674 


4.93158 


5.02609 


6 


5.44414 


5.50073 


5.91700 


5.99970 


7 


6.44353 


6.48362 


6.91468 


6.97787 


8 


7.44614 


7.47397 


7.91743 


7.96350 


9 


8.44969 


8.46922 


8.92199 


8.95537 


10 


9.45331 


9.46734 


9.92696 


9.95146 


11 


10.4567 


10.4670 


10.9318 


10.9501 


12 


11.4598 


11.4676 


11.9362 


11.9503 


13 


12.4625 


12.4685 


12.9403 


12.9512 


14 


13.4649 


13.4697 


13.9440 


13.9527 


15 


14.4671 


14.4709 


14.9472 


14.9543 


16 


15.4690 


15.4721 


15.9502 


15.9560 


17 


16.4707 


16.4733 


16.9529 


16.9576 


18 


17.4723 


17.4745 


17.9553 


17.9593 


19 


18.4737 


18.4755 


18.9575 


18.9609 


20 


19.4750 


19.4765 


19.9595 


19.9624 


21 


20.4762 


20.4775 


20.9613 


20.9638 


22 


21.4772 


21.4784 


21.9630 


21.9651 


23 


22.4782 


22.4792 


22.9645 


22.9664 


24 


23.4791 


23.4800 


23.9659 


23.9676 


25 


24.4799 


24.4807 


24.9672 


24.9687 


26 


25.4807 


25.4814 


25.9685 


25.9698 


27 


26.4814 


26.4820 


26.9696 


26.9708 


28 


27.4820 


27.4826 


27.9706 


27.9717 


29 


28.4827 


28.4831 


28.9716 


28.9726 


30 


29.4832 


29.4837 


29.9726 


29.9734 


31 


30.4838 


30.4842 


30.9734 


30.9742 


32 


31.4843 


31.4846 


31.9742 


31.9749 



Table 2: The values of the improved distances r\/a Eq. (lA.ip and r/a Eq. (|A.5j) extrap- 
olated to L/a — >• oo for the case of HYP and HYP2 smearings. We show 6 significant 
digits for all values of r/a, where the last digit is rounded. 
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(ism = 1 for unsmeared static quark lines). 

The distance f defined from Eq. (12.171) is given in the case of HYP smeared 
static quarks by 

r 3 = 2Tr[GuYp(r + a) + G HYP (r-a)-2GuYp{r)}/a 2 . (A.5) 

In practice we evaluate the momentum integral in Eq. ( 1A.2I) by discrete momentum 
sums over pj = 2nnja/L, rij — 0, 1, . . . , L/a — 1. In Table |2] we quote the results 
for T\ja and f/a obtained from extrapolations L/a —t oo. The latter are done with 
the method explained in Appendix D of [51] and we consider lattice sizes larger 
than L/a = 128 up to L/a = 512. Due to the symmetry under pj — > —pj of the 
integrand only odd powers of a/L can appear in the expansion in powers of a/L 
and in general this evaluation of the integral is the application of a trapezoidal 
rule, which has discretization errors of order (a/L) 2 . Thus the leading correction 
is S\(a/L) 3 . The data for r\ja and f/a are very well fitted by a polynomial 
so + Si(a/L) 3 + S2(a/L) 5 and we added terms S3 (a/L) 7 + s^a/L) 9 to estimate the 
systematic error of the extrapolations. In Table |2] we list the extrapolated values 
with six significant digits. 



B Perturbation theory for c(r) 



We consider QCD with iVf massless dynamical quark flavors. The quantity c(r) 
in Eq. (I1.2p defines a renormalized coupling (Cf = 4/3), 

4:71 



9 C 0) 



-c(r) , pL — l/r. 



(b.i; 



It is very similar to (7q q (/i) = ^r 2 F{r) , /i = 1/r discussed in [50]. The relation is 

For a perturbative evaluation of one-scale quantities such as c(r) it is natural to 
just integrate the renormalization group equation 



d/i 



(B.3) 



We do this in the precise form of 



A, 



-2 r 6i/(2bg) e -l/(26 ^) exp 



dx 



(3 c (x) b x 3 b\x 



1 1 

+ 



(B.4) 



8 It has also been observed in more than one case that it also yields a good perturbative 
description of the non-perturbative behavior. 
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where for f3 c the truncated perturbative expansion is inserted, but the integral is 
(numerically) evaluated as it stands. The Lambda-parameter in the c-scheme is 

A c = e- 1 /2 Aqq = A_ efcl /(8.M-i/2 ) (R5) 

where fa = ^(ai + a 2 N f ), ai = -35/3 + 22 7e and a 2 = 2/9 - 4 7e /3 [55J . We 
now turn to the perturbative beta function. 



B.l Perturbative beta function in the c-scheme 

The expansion of the potential in the MS coupling is now known to a high accuracy. 
After the term [56j[57], the resummation of the infrared divergent diagrams 
appearing first at the next order was performed [521158] . yielding a ~ g|^log(g|jg-) 
term. Recently also the full three-loop computation was finished by two groups 
[SHIED]- Due to the log (g^) term in the potential [521I2E], the beta function 
has a perturbative expansion 

= -f c iit b n9c n + bf}f c \og(C A f c /(Sn)) + O(^)] , (B.6) 

n=0 

with the universal coefficients (C\ = 3) 

b^ = b = -i-(llC A /3-2iVf/3), (B.7) 



(4tt) 2 ' 

bf=h = ^(34^/3-10^^/3-2^^). (B.8) 

We now describe how the non-universal coefficients are obtained from the 
results in the literature. Our starting point is Eq. (40) of [51], which is the 
expansion of the static potential V(r), denoted "static energy" in [5T], in the 
MS strong coupling a s = g^(l / r) / (An) derived from the above mentioned work. 
Introducing the notation V(r) = — CpG(a s )/r, we obtain an expansion for a c = 

a c = l -r 2 G'\a s )-rG\a s ) + G(a s ) (B.9) 

= a s + di a 2 s + d 2 0% + d 3 a 4 s + d 3>l « s 4 In (^y^J + 0{a\) , (B.10) 

where the primes in the first equation mean derivatives with respect to r and the 
coefficients of the expansion are 

di = -jj- (Si - 36 (4tt) 2 ) , (B.ll) 



20 



d 2 
d 3 

d 3 ,i 



(An) 2 
1 



(a 2)S + 46 2 (47r) 4 - 36!(47r) 4 - 65 1 6 (47r) 2 ) 
(a 3 + 12a 1 6g(47r) 4 - 6hh (4tt) 4 + 106 &i(47r) f 



(4tt) 3 

-36 2 (47r) 6 -9a 2iS 6 (47T) 2 ) 



12tt 



(B.12) 

(B.13) 
(B.14) 



The coefficients b 2 and b 3 of the beta function in the MS scheme, /3ms(^ms) ~ 
— ^Ms"Sn>o ^m§"^™' can ^ e f° un d m [621163] - The coefficient di is defined in Eq. (7) 
and the coefficient a 2jS in Eq. (8) of [61], they both depend on Nf. The coefficient 
a 3 is 

«3 = 7 ^(c (iV f ) + 27 E c 1 (iV f ) + (47 2 + 7r 2 /3)c 2 (iV f ) 



3CY 



+ (8 7 I + 27r 2 7 E + 16C(3))c 3 (iV f )) 



where 



c (N { 



219.59 + ( 4 1} iV f + a^'N'f + a^iV 3 ) /4 



(2) 



,(3) 



(B.15) 



(B.16) 



From [61JE1] we getfH Co(0) and the coefficients a^, a^ 1 and a® are given in Eq. 
(6) of [65]. The coefficients ci(iVf), c 2 (iVf) and c 3 (Nf) are defined in Eqs. (10), 
(11) and (12) of [66] respectively. 

The non-universal coefficients b^ and b^ as well as the coefficient b^\ may 
now be computed by differentiating (3 C = y^^p with a c of Eq. (IB.lOp . where 
the MS beta function is used. This first yields /3 C as a function of a s from 
which we change to (3 c (g c ) by inserting the inverted Eq. (IB.lOj) . a s = a c + . . . — 
4,«« c 4 ln(^). 

Carrying this out in MAPLE we find 



.(c) 



<(c) 



b 2 - 56 3 + fi 2iS 6 (4vr)- 4 - &A(47r)- 2 - a 2 6 (4vr)- 4 
(4tt)" 3 [0.98165 - 0.16738iV f - 0.002l2iV f 2 + 0.000267V f 3 ] 



6 3 - 2a 1 6 2 (4 



7T 



i-2 



J 3,l 



+ 2a 3 & (4vr)- 6 + -C 3 6 (47r)- 4 - 256g&! 

-6& 1 a 2iS 6 (47r)- 6 + dlhiAn)- 4 - 36b 4 + 4^6 (4tt)" 6 

(4tt)- 4 [0. 12206 + 0.09696iV f - 0.01899iV f 2 + 0.0004458iV f 3 

+0.0000195iV f 4 ] 
2 

-Ci& (47r)~ 4 = (4vr)- 4 [1.25385 - 0.07599iV f ] . 



(B.17) 



(B.18) 



(B.19) 



3 We thank the authors of [61] for communication on the value of cq(0). 
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Figure 8: The perturbative running of the quantity c(r) obtained from Eq. (|B.4p using 
the 2-loop (dotted lines), 3-loop (dashed lines) and 4- loop (continued lines) beta function 
/3 c (5c) for iV f = 0,2. 

As in the MS scheme the coefficients (4ir) n+1 b n are of order one and thus the beta 
fuction has a well-behaved expansion up to couplings a c of the order of 1/3. The 
perturbative running of c(r) is shown in Fig. [HJ 

The "asymptotic convergence" of the series Eq. (IB.lOp is not good. It can 
be substantially improved by matching the couplings at a different scale, i.e., 
by expressing a c (s/r) as a function of a s (l/r) and choosing s = so = A c /Aj^ 
("fastest apparent convergence", cf. [501 [55]). The resulting curves for a c are 
hardly distinguishable from the ones shown in Fig. [8j 

B.2 Perturbative beta function in the qq-scheme 

In the same way one obtains the beta function in the qq-scheme. We update the 
formulae given in [55] to include the 4-loop term: 




4 qq) = (47r)" 3 [1.6524- 0.28933iV f + 0.00527N f 2 + 0.00011iV f 3 ] (B.20) 
6^ qq) = (4tt)- 4 [4.94522 - 1.07965iV f + 0.079107^ - 0.002774iV f 3 




+0.000051iV f 4 ] 



(B.21) 
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2 

b^f = -Ci&o(47r)~ 4 = (4tt)- 4 [1.25385 - 0.07599^ f ] . (B.22) 

Perturbation theory in the c-scheme appears much better behaved than in the 
qq-scheme. Since this can only be considered an accident we come to the same 
conclusion as [22], namely that the perturbative description of the static poten- 
tial is accurately valid only at rather small values of r, where a qq (l/r) ~ 1/4. 
Unfortunately these distances are close to present lattice spacings. In particular 
the data presented in this paper are not good enough to extract the A parameter 
through Eq. ( 1B.4j) or variants thereof. 
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